Cytoscape可视化物种分类树结构

前两天偶然看到了这样一张图，听群友说这是Cytoscape调出来的效果，就想试一试。

折腾一番后，发现这个图的工作量貌似蛮大的。因此本篇就只是概括了一下大致的过程，包括如何准备输入数据，以及在Cytoscape可视化时大致上要点击哪些按钮等，并未仔细琢磨成品图。所以，下文略显粗糙，大家有兴趣就看吧……

这其实是个物种分类树，代表的微生物分类。从高级分类水平（界，细菌和真菌）的根节点开始，依次往低分类层级（门纲目科属等）延伸，逐渐形成分支，末端节点应该是OTU或者种水平。因此，分支呈簇状展示，分支的颜色主要按门水平着色，将衍生自同一类门的微生物分支赋值相同的颜色。节点的大小应该对应了微生物的丰度，而节点的颜色主要有两种，猜测红色代表了一些重要的微生物，白色是不重要的。图中文字部分则标识了主要的微生物门所含的节点数量，以及这些微生物的总丰度占比。

好吧，其实白鱼同学并不知道这是来自哪篇文章的插图，不清楚它的任何细节，包括其具体的指代信息以及作者试图表达的意义等，上一段文字纯粹是看图瞎猜的……又懒得问其来源，那就按照这种理解去模仿吧，反正只是随手画个图嘛，不要在意那么多细节。

准备数据

示例数据和预处理R代码的百度盘链接（提取码，sf9f）：

https://pan.baidu.com/s/1nS402FLH3lQu8CNq6JjEsA

网盘附件“otu.txt”如下所示，以OTU的丰度为主，后面跟上这些OTU所属的界门纲目科属分类。需要事先把一些注释不明确的物种提前去除（比如unidentified、unknown、others等），否则导致了一些分类关系混乱。然后将它读到R中，整合各个分类水平的包含关系，即上一级分类和下一级分类的对应列表用作边列表，并统计共计多少分类类型及其丰度用作节点列表。

#读入数据
otu <- read.delim('otu.txt', stringsAsFactors = FALSE, check.names = FALSE)

#整合各分类层级和所含物种丰度的嵌套关系，构建网络边列表
#source 和 target 分别对应了上下级物种分类层级的对应关系，abundance 表示 target 的丰度
genus_otu <- otu[c('genus', 'otu', 'abundance')]
names(genus_otu) <- c('source', 'target', 'abundance')
family_genus <- aggregate(otu$abundance, by = list(otu$family, otu$genus), FUN = sum)
names(family_genus) <- c('source', 'target', 'abundance')
order_family <- aggregate(otu$abundance, by = list(otu$order, otu$family), FUN = sum)
names(order_family) <- c('source', 'target', 'abundance')
class_order <- aggregate(otu$abundance, by = list(otu$class, otu$order), FUN = sum)
names(class_order) <- c('source', 'target', 'abundance')
phylum_class <- aggregate(otu$abundance, by = list(otu$phylum, otu$class), FUN = sum)
names(phylum_class) <- c('source', 'target', 'abundance')
kingdom_phylum <- aggregate(otu$abundance, by = list(otu$kingdom, otu$phylum), FUN = sum)
names(kingdom_phylum) <- c('source', 'target', 'abundance')

edge_list <- rbind(kingdom_phylum, phylum_class, class_order, order_family, family_genus, genus_otu)

#将各个分类层级对应至 phylum，以便用于在 Cytoscape 中调整树分支的颜色
otu_phylum <- otu[c('phylum', 'otu')]
names(otu_phylum)[2] <- 'target'
genus_phylum <- otu[c('phylum', 'genus')]
names(genus_phylum)[2] <- 'target'
family_phylum <- otu[c('phylum', 'family')]
names(family_phylum)[2] <- 'target'
order_phylum <- otu[c('phylum', 'order')]
names(order_phylum)[2] <- 'target'
class_phylum <- otu[c('phylum', 'class')]
names(class_phylum)[2] <- 'target'
phylum_phylum <- otu[c('phylum', 'phylum')]
names(phylum_phylum)[2] <- 'target'

tax_phylum <- rbind(phylum_phylum, class_phylum, order_phylum, family_phylum, genus_phylum, otu_phylum)
tax_phylum <- tax_phylum[!duplicated(tax_phylum$target), ]
edge_list <- merge(edge_list, tax_phylum, by = 'target', all.x = TRUE)

#输出边列表
edge_list <- edge_list[c('source', 'target', 'abundance', 'phylum')]
edge_list$weight <- 1
head(edge_list)

write.table(edge_list, 'edge_list.txt', row.names = FALSE, sep = '\t', quote = FALSE)
#构建节点列表，包括各个分类水平的名称，以及所含物种总丰度等
node_list <- reshape2::melt(otu, id = 'abundance')
node_list <- aggregate(node_list$abundance, by = list(node_list$value, node_list$variable), FUN = sum)
names(node_list) <- c('shared name', 'taxonomy', 'abundance')

#原始丰度 abundance 的数值差较大，不利于直接使用该列在 Cytoscape 中定义节点大小
#再指定个新列 abundance2，对原始的 abundance 作些转换，比如 log 降低高丰度的权重等
node_list$abundance2 <- log(node_list$abundance+1, 2) #log 转化时加个基数 1

#添加节点的 phylum 水平分类
otu_phylum <- otu[c('phylum', 'otu')]
names(otu_phylum)[2] <- 'shared name'
genus_phylum <- otu[c('phylum', 'genus')]
names(genus_phylum)[2] <- 'shared name'
family_phylum <- otu[c('phylum', 'family')]
names(family_phylum)[2] <- 'shared name'
order_phylum <- otu[c('phylum', 'order')]
names(order_phylum)[2] <- 'shared name'
class_phylum <- otu[c('phylum', 'class')]
names(class_phylum)[2] <- 'shared name'
phylum_phylum <- otu[c('phylum', 'phylum')]
names(phylum_phylum)[2] <- 'shared name'

tax_phylum <- rbind(phylum_phylum, class_phylum, order_phylum, family_phylum, genus_phylum, otu_phylum)
tax_phylum <- tax_phylum[!duplicated(tax_phylum$'shared name'), ]
node_list <- merge(node_list, tax_phylum, by = 'shared name', all.x = TRUE)

#根据事先选择的一些重要节点名称，在节点列表中标识出后输出节点列表
select_node <- read.delim('select_node.txt', stringsAsFactors = FALSE)

node_list[which(node_list$'shared name' %in% select_node$id),'select'] <- 1
node_list[which(! node_list$'shared name' %in% select_node$id),'select'] <- 0
head(node_list)

write.table(node_list, 'node_list.txt', row.names = FALSE, sep = '\t', quote = FALSE)

Cytoscape的可视化

这样两个输入文件就准备好了，分别为上述R输出的“edge_list.txt”和“node_list.txt”，导入至Cytoscape中进行可视化。

1 读取输入文件

打开Cytoscape后，点击“Import Network from File System”读取边列表，也就是上述输出的“edge_list.txt”；以及点击“Import Table from File”读取节点列表，“node_list.txt”。

文件导入后，Cytoscape中自动呈现一幅网络图，这个网络图的结构就代表了给定的物种分类树了。

2 边的颜色和尺寸等属性调整

点击界面左侧“Style”按钮进入外观调整选项。首先是对边的调整，点击下方“Edge”后进入边样式调整界面。

已知该网络代表了物种分类树，如上所述，从高级分类水平（界）的根节点开始，依次往低分类层级（门纲目科属等）延伸，逐渐形成分支。因此，下游的分类分支都衍生自上游高级分类分支。参考图中将衍生自同一门分类的微生物分支赋值相同的颜色，这里我们效仿，在上文的R操作过程中已经对网络中每条边做了归类，边列表中的“phylum”列就记录了分支对应的门分类水平。点击“Stroke Color”，按边列表中的“phylum”列分配颜色。

边的颜色是主要的调节属性。对于其它的属性，视情况自定义修改。例如如果觉得边太细，就可以点击“Width”将边的尺寸设置的宽一些，使其更清晰。

3 背景色调整

随后，点击下方“Network”，将“Background Paint”背景色设置为黑色。

4 节点的颜色和尺寸等属性调整

然后是对节点的调整，点击下方“Node”后进入节点样式调整界面。

对于节点颜色，点击“Fill Color”，按节点列表中的“select”列分配颜色，将1（标记的重要的微生物）赋值为红色，0（不重要的）赋值为白色。

对于节点形状，勾选下方“Lock node width and height”后，点击“Shape”统一设置为圆形。

对于节点大小，点击“Size”，按节点列表中的“abundance2”列指定大小。该列是个连续的数值，大值对应大点，小值对应小点。

以及点击“Properties”显示出“Label Transparency”选项后统一设置为0，即将节点标签均设置为透明，以去除节点标签。

5 网络布局调整

上述过程都很简单，不再多说了，大家看着来设置就行。

最后是调整网络布局，获得类似示例图中的树状结构。可以先尝试一些自动布局，看看哪些布局风格效果较好。下图展示一些和示例图较为相似的自动布局，以及一些不相似但是可能有趣的自动布局，供大家参考吧。

点击菜单栏“Layout”，提供了N种预设布局；此外“Apps”中可以下载更多样式的插件，例如yFiles布局就是在“Apps”中找的拓展插件之一。

这些自动布局也就是有几分相似而已。一般来说，网络图布局光靠自动布局是难以满足需求的，配合手工调整肯定是免不了的，只是手工调整的过程会很糟心……哦对了，说不定真的有自动布局能够出来示例图样式，只是没找到而已。

以下简单展示一个手工调节布局的过程。以“Prefuse Force Directed Layout （按 weight 列定义布局）”布局为例，因为它的微生物节点分类层次不明确，有待手动拖一下。当然，实际情况中尽可能找一个最相似的布局来，这样可以简化手动的工作量。点击“Select”，按微生物节点所属的门分类水平进行选择，并将选中的节点拉动至一旁，就可和其它类群分开了。随后，设法手动选择节点子集，并拖动将该类群内部的节点的层次结构尽可能区分明显。

把所有节点的层次关系调整好，将微生物的分类关系较好地呈现出，总之工作量挺大的。所以，白鱼同学就没再继续，为了一个示例耗费大量时间在里面就不值得了……本篇教程只是作为方法类指引，帮助大家了解大体上是怎样的操作过程就可以了。真正有需要时，再花点时间配合手工琢磨布局，再慢1天也能搞定了吧。

6 关于文字图例的添加

最后，如果在Cytoscape中调试出了想要的结果，点击“File > Export > Network to Image”将图片导出为pdf矢量图后，放在AI（Adobe Illustrator）中补充文字标签。例如示例图中，微生物门水平的名称、节点数量和所含物种的丰度占比等，就可以通过AI补充。

整体过程大致就这样吧，最麻烦的在布局调整那一步，其它倒还好。总之真要有待于可视化这样一张图的话，有耐心就是了……集中精力，1天差不多也能搞定了。

往期精品(点击图片直达文字对应教程)

后台回复“生信宝典福利第一波”或点击阅读原文获取教程合集